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CIh, Type la supernovae are thought to be caused by thermonuclear explosions of a 

Q I carbon-oxygen white dwarf in close binary systems. In the single-degenerate scenario 

-t^ ' (SDS), the companion star is non-degenerate and can be significantly affected by the 

C^ ■ explosion. We explore this interaction by means of multi-dimensional adaptive mesh 

refinement simulations using the FLASH code. We consider several different companion 

CN , types, including main-sequence-like stars (MS), red giants (RG), and helium stars (He). 

^sq , In addition, we include the symmetry-breaking effects of orbital motion, rotation of the 

cn I non-degenerate star, and Roche-lobe overflow. A detailed study of a sub-grid model 

for Type la supernovae is also presented. We find that the dependence of the unbound 

^f-N ■ stellar mass and kick velocity on the initial binary separation can be fitted by power-law 

^D ' relations. By using the tracer particles in FLASH, the process leading to the unbinding 

of matter is dominated by ablation, which has usually been neglected in past analytical 

studies. The level of Ni/Fe contamination of the companion that results from the 

L^ , passage of the supernova ejecta is found to be a '^ 10~^Mq for the MS star, ~ 10~^Mq 

^ [ for the He star, and ~ 10~®Mq for the RG. The spinning MS companion star loses 

- . .' about half of its initial angular momentum during the impact, causing the rotational 

velocity to drop to a quarter of the original rotational velocity, suggesting that the 

Tycho G star is a promising progenitor candidate in the SDS. 

Subject headings: binaries: close — supernovae: general — methods: numerical 



1. INTRODUCTION 

Type la supernovae (SNe la) are consid ered to be "standardizable" candles when measuring 



the distance of distant objects (jPhillipall993l ) and thus play an important role in cosmology, helping 



1998; Perlmutter et al 



determine key cosmological p arameters (JBranch &: Tammannlll992l : lRiess et al.lll998l : ISchmidt et al 



19991 ) . Most cosmological applications using SNe la are based on an assump- 



tion that the homogeneity of SNe la applies to high-redshift SNe la. However, this assumption is 
questionable until we fully understand the nature of the progenitors of SNe la. According to recent 
studies, SNe la are believed to be thermonuclear explosions of carbon-oxygen (CO) white dwarfs 
(WDs), but the ignition mechanism and the progenitor systems are still unknown. 

Current p opular progenitor mode ls are class i fied in to two major scenarios, the single-degenerate 
scenario fSDS.IWhelan k. IbenI (119731): iNomotd (jl982l )) and the double-degenerate scenario (DDS, 
Iben fc Tutukovl (|l984l ): Iwebbinki (|l984l )). The SDS involves a CO WD accreting matter from a 
non-degenerate binary companion, such as main-sequence (MS) stars, red giants (RG), or helium 
stars (He). The CO WD eventually becomes unstable and then, when the mass reaches the Chan- 
drasekhar limit (M ~ IAMq), it explodes as a SN la. The DDS instead involves two CO WDs, 
with total mass larger than the Chandrasekhar mass, that merge together because of the emission 
of gravitational waves. 

With the DDS, the lack of hydrogen spectral lines is easy to explain; however, the homogeneity 
of the explosion energy is difficult to account for since the total mass of the binary system varies 
from case to case. Furthermore, it is possible that WD- WD mergers end with the production o f 
an ONeMg WD followe d by accretion- induced collapse to a neutron star (jNomoto &: IbenI Il985l ). 
Badenes &: Maoa (|2012l ) calculate the merger rate of binary WDs in the Galactic disk based on 
the observational data in the Sloan Digital Sky Survey. They conclude that the merger rate of 
binary WDs with super-Chandrasekhar masses would not significantly contribute to the SNe la 
rate. The SDS produces a relatively fixed explosion energy, which accounts for the homogeneity. 
The non-degenerate companions in the SDS, however, are usual ly hydrogen-rich which leads to the 
H/He contamination problem in the SDS (see iLivid ()200d ) and iHillebrandt fc Niemeverl (J200d ) for 
detailed descriptions of the advantages and weaknesses of the SDS and DDS). 

Recent studies based on the delay time dist ribution (DTD) suggest the ne e d for at least a 



two-co mpo nent model for the D TD. In particular, IScannapieco fc Bildstenl ()2005l ). lMannucci et al. 



(|200d ) and lBrandt et al.l (l20ld ) found that the observations can be fitted with a short-delay-time 
population (~ 10® yr) and a long-delay-time population (3 — 4 Gyr). 

The calculation s of SN la rates and DTD using binary populat ion synthesis c an be matched 



2008 



with both the SDS (Han fc Podsiadlowski 



2003 



Hachisu et al.ll2008l ) and the DDS (|Pritchet et al.l 



Ruiter et al.ll2009l ). The uncertainties in the observed DTD, howe ver, are dominated by th e 



uncertainties in galactic stellar populations and star formation histories ( Maoz fc Mannuccill201ll ). 
Thus, it is still hard to distinguish between the viability of the SDS and the DDS progenitor models 
without more detailed observations. 

With detailed binary calculations, the long-delay-time population can be und erstood in terms 
of progenitor systems characterized by a M S-like companion in the MS- WD channel (llyanova &: Taam 



2004 



_Dy ; 
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Wang et al.ll20ld : IWang fc Hanll20ld ) and/or by an RG in the RG-WD channel (JHachisu et al. 



1999I . |2008| ). In contrast, the shor t-delay-time popula tion may consist of systems with a massive 
MS star in the M S- WD channel ([Hachisu et al.l l2008l ) or with a He-star in the He- WD channel 
(|Wang et al.1 12003 ) . Therefore, the more likely explanation is that both SDS and DDS have con- 
tributed to SNe la. 

Several detailed h ydrodynamics simulations with the SDS have been studied in the past decade. 
Marietta et al.l (J200d ) examined the impac t of a SN la on a MS star, a sub- giant, and a RG using 
two-dimensional Eulerian hydrodynamics. [Marietta et al.l (J200d ) found that the MS star and sub- 
giant companions lose 15% of their mass after the explosion, and the RG companion loses about 
96% — 98% of its envelope. They also found that the impact of SN ejecta with the companion star 
creates a hole with an opening angle of ~ 30° in the high-velocity ejecta and ~ 40° in the low- 
velocity ejecta. This hole corresponds to 7% — 12% of the ejecta's surface, making the supernova 
remnant (SNR) anisotropic and potentially affecting Si H spectral line shapes. 



Pakmor et al.l (|2008l ) re-examined the M S- WD channel arid upd ated these results by consider- 



ing the pre-supernova binary evolution from llvanova &: TaamI (J2004l ) and using three-dimensional 
smoothed particle hydrodynamics (SPH) simulations. They found strong dependences of unbound 
mass and kick velocity on initial binary separation with power-law relations. Because of the rel- 
atively small number of part i cles i n their SPH simulations, they were unable to reprod uce the 
turbulence in [Marietta et al.l (120001 ). However, their comparison to [Marietta et al.l (|2000l ) found 
that the unbound mass does not depend much on the numerical method used. Although the 
unbound mass is not sensitive to the turbulence around the MS star, other important physical 
quantities, e.g. contaminated SN ejecta on the companion star, may be sensitive to small-scale 
turbulence. 



In cont r ast w ith [Marietta et al.l (J2000l ) , |Pan et al.l (|20ld ) studied the He- WD channel from 
Wang et al.l ()2009l ) for the short-delay-time population, using two-dimensional Eulerian hydrody- 
namics simulations with adaptive mes h refinement. They a lso found the same power-law relations 
in unbound mass and kick velocity as iPakmor et al.l (|2008l ) . but with different power-law indices. 
An upper limit of the ratio of Ni/Fe contamination to the helium abundance was found to be 
~ 9 — 50 X 10"'^, which is higher than the solar abundance ratio. This result suggests the possibility 
of detecting remnant helium stars in their SNR. 

In this paper, we examine the SDS using the above discussed companion star models to de- 
termine whether the contamination problem of H/He in the companion envelope can be overcome. 
We revisit this problem using modern adaptive mesh refinement (AMR) and different companion 
models that incorporate the effects of binary evolution. We include the symmetry-breaking effects 
of orbital motion, rotation of the non-degenerate star, and Roche-lobe overflow (RLOF). In the 
next section, the numerical methods and the construction of initial setups are described. Our nu- 
merical results for different companion models are reported in Section 3. In Section 4 we present 
a parameter study made by varying significant physical quantities and then discuss the possible 
observational implications. In the final section, we summarize our results and conclude. 



NUMERICAL METHODS AND MODELS 



2.1. Numerical Codes 



Two different codes are employed in this work. The first one is a stellar evolution code used to 
construct the non-degenerate stellar models. The second code is a multi-dimensional hydrodynamics 
code to simulate the impact of type la supernova ejecta. 

The stellar evo lution code used is MESAlU (Modules for Experiments in Stellar Astrophysics, 
Paxton et al.l (J201ll )). MESA combines many sophisticated numerical methods to control time and 
space resolutions, and many state-of-the-art physical modules to solve the stellar structure. For a 
given initial stellar mass and metallicity, MESA evolves stars to a certain age and provides all the 
physical quantities we need for the hydrodynamics simulations. 

For the SN la simulations, we use FLASfi^l version 3 (JFryxeU et ahlboool : bubey et allboosl ). 
FLASH is a parallel, multi-dimensional hydrodynamics code based on block-structured adaptive 
mesh refineme nt (AMR). To solve the Eul er equations, we use the split piecewise parabolic method 



Helmholtz EOS (JTimmes fc Swestv 



(PPM) solver (ICoMla fc Woodward Il984l ) in FLASH. The equation of state (EOS) applied is the 



20001 ) ■ which is interpolated from a precomputed table of the 



Helmholtz free energy. It includes contributions from radiation, completely ionized nuclei, and 
degenerate electrons and positrons for an optically thick mixture of gas and radiation in local 
thermodynamic equilibrium. Magnetic fields are ignored, but self-gravity is considered and solved 
using the multipole Poisson solver in 2D and using the multigrid Poisson Solver in 3D (JRicker 
20081 ). Particle modules in FLASH 3 have the ability to treat both active and passive particles. 
Active particles are massive particles that interact through gravity with other active particles and 
fluid elements. Passive particles are massless particles that only follow the motion of fluid elements 
in Lagrangian coordinates. We use active particles to represent the CO WD and the core of the 
RG companion and passive particles to monitor the motion of fluids. 



2.2. Non-Degenerate Companion Models 



In order to compare our setup with previous work by [Marietta et al.l ()2000l ) and iPakmor et al. 



(120081). we use MESA to create a MS star c ompanion mode l simi lar to the HCV scenario in 
Marietta et al.l (120001 ). The HCV scenario in [Marietta et al.l (J2000l ) is a hydrogen cataclysmic 
variable system consisting of a MS star and a CO WD. We assume the binary system is in RLOF 
and has a binary separation a = 3i?* at the onset of the SN la explosion. The MS star companion 
in the HCV scenario has a mass M = 1.017 Mrr, and radius i?* = 6.8 x 10^ km. Because of the lack 



^http : //mesa . sourcef orge . net/index . html | 
^http : //flash . uchicago . edu| 
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Fig. 1. — Mass vs. radius profiles for the main sequence, red giant, and helium star models. For 
the RG-WD scenario, only the envelope mass of the RG is plotted here. 



of certain physical data in the HCV scenario, we created a similar MS star model, named HCV, 
with M = 1.027 Mq, and R^ = 6.76 x lO^km. 

For more realistic models, we use models of SN la progenitor candidates taken from detailed 
binary calculations. Based on current studies, we chose three possible stellar types of binary 
companions in the SDS: MS stars, RGs, and He stars. 



Ivanova &: TaamI (J2004l ) studied the evolution of binaries consisting of evolved MS stars with 



WDs. They also investigated possible channels that may eventually evolve to SNe la. We chose 
parameters from the case with final mass Afd;f = 1.17M0 for the MS star and initial orbital period 
P = 1 day to construct a model of a MS- WD binary system (the system is denoted by MS-WD in 
our simulations). 



Hachisu et al.l (jl999l ) have proposed an alternative binary channel consisting of a WD and a 



low-mass RG. The accreting WD in this model has a strong optically thi ck wind, which broad ens 
the range of parameters leading to SNe la. We also adopt the case PI in lHachisu et al.l (jl999l ) for 
our simulations (denoted by RG-WD here). This system consists of a one solar mass WD and a 
two solar mass RG with an initial orbital period Pq = 300 days. The system ends with a SN la 
after 7.2 x 10^ years; the mass of RG becomes 0.98Mq at the onset of the SN la explosion. 



For the helium star model, we rely on the model originally proposed by I Wang et al.l (J2009l ) 
as well as the case where t he binary system is in a stable He-shell burning phase at the onset of 
the SN explosion (Case 2 in lWang et al.l (J2009l )). This model also is the He-WDc model applied in 
Pan et al.l (120101 ) . but with nonuniform chemical composition. The mass of this helium star model 
is 1.OO7M0 (denoted by He-WDc). 

In MESA, we do not follow full binary evolution but use the initial conditions with a constant 
mass-loss rate, estimated from the above references as an approximation. The mass-loss rate used 
is taken to be the ratio of the mass change to the evolution time in the above references. All models 
we created from MESA are summarized in Table [TJ and the mass versus radius profiles are shown 
in Figured) 



2.3. Initial Setup 

We conduct two- and three-dimensional hydrodynamics simulations to investigate the impact of 
SN la ejecta on the binary companions. In the three-dimensional calculations, we rely on Cartesian 
coordinates; the simulation box dimensions are set to a size equal to 30 times the radius of the non- 
degenerate companion star {R^,) in all directions. The companion star is located at the center of the 
simulation box for convenience. Two-dimensional models are calculated in cylindrical coordinates 
due to the axial symmetry if the orbital motion of the binary system is ignored. The simulation 
box dimensions are set to IbR.^ in the radial (r) direction and 30-R,, in the axial [z) direction. We 
interpolated the one-dimensional model onto the FLASH grids using up to 10 levels of refinement 



based on the magnitudes of the second derivatives of gas density and pressure. Each AMR block 
contains 8^ zones in the three-dimensional box and 8 x 16 zones in the two-dimensional box. This 
corresponds to an effective uniform resolution of 4096^ in three-dimension and 4096 x 8192 in two 
dimensions. To simplify the problem, the compositions from MESA are adjusted to hydrogen (^H), 
helium (^He) and carbon (^^C) only (companion material). We use outflow boundary conditions 
for fluids and isolated boundary conditions for the Poisson solver. 

In an Eulerian hydrodynamics code like FLASH, physical properties are calculated at fixed 
spatial positions. It is difficult to trace fluid elements in a time sequence. However, FLASH has the 
ability to trace fluid elements using passive particles. Passive particles are massless particles that 
only follow the motion of fluid elements in Lagrangian coordinates without interacting with fluids. 
In our simulations, 10^ passive particles are distributed with the gas density in order to study the 
history of shocked gas. 

For the RG model, although we have AMR, the core of the RG still cannot be resolved with 
current computational resources. Thus, we use a rigid spherical particle cloud with 3 x 10^ active 
particles to represent it. We artificially modify the physical quantities within the core region of the 
RG to increase linearly with the radius and then reconstruct the model on multi-dimensional grids 
in FLASH. Therefore, the mass difference between the real stellar model and the reconstructed 
model in FLASH is the core mass of the RG (Mfj,G,c = O.SIIMq). The radius of the particle cloud 
is set to three times the smallest zone spacing; the force on the cloud uses cloud in cell (GIG) 
interpolation. We utilize a particle cloud , instead of a single particle, to avoid the problem of force 



anisotropy in GIG (JRicker &: Taainll2008l ). 



In order to validate our models in hydrostatic equilibrium and to reduce the geometrical 
distortion between one dimension and multi-dimensions, all models are relaxed on the multi- 
dimensional Eulerian grid by artificially damping the gas velocity for 5 dynamical time scales 
(^dyn = l/2(Gp)~^"). During this time, a damping factor, smoothly increased from 0.7 to 0.99 
on the gas velocity, is imposed at each timestep. We ensure that the companion models remain 
relatively stable during the relaxation and that the Mach number is always smaller than 0.01. 

After relaxation, we introduce asymmetric effects from orbital motion by adding a WD into 
the simulation box (three-dimensional runs only) . The WD is represented by another particle cloud 
with 3 X 10^ particles and its mass is set to M = 1.378M0 (we assume the binary system is about 
to explode). In the three-dimensional simulations, the orbital plane is set on the x — y plane and 
the WD is placed on the positive x— axis (positive radial direction in two dimension) with a binary 
separation a. Positive z— axis is set to the direction of angular momentum. The binary system is 



assumed to be in RLOF and the binary separation is calculated using equation [T] (lEggletonlll983l ): 

Rl 0.49g2/3 



a 0.6g2/3 + ln(l + gV3) ' 



(1) 



where Rl is the radius of the Roche lobe-filling star, and q is the mass ratio. The non-degenerate 
binary companion has a spin with a spin-to-orbit ratio of 0.95. In this phase, the maximum AMR 



level is reduced to 8 levels to save computation time and the AMR level is forced to the maximum 
level (8) in two spherical regions centered on the WD and the companion star within a radius of 
R < l.SR^,. AMR levels in other regions are calculated from the second derivative of gas density 
and pressure and limited to 6 levels (a situation denoted by 6/8 levels, equivalent to a 1024^ uniform 
grid in three dimensions and 1024 x 2048 in two dimensions). 



2.4. The Supernova Model 

After the binary system is placed in the simulation box, we simulate the binary evolution for a 
few orbital periods. During this phase, mass transfer via RLOF occurs and the companion star is 
slightly distorted due to the tidal force. With current computational resources, we cannot simulate 
the whole delay time in multi-dimensional hydrodynamics simulations. We therefore simulate the 
binary evolution for only a few orbital periods (< 3Porb) and assume the morphology does not 
change much once the RLOF is stable. We tested the impact of a SN la explosion on the binary 
companion using different binary evolution times. We find that binary evolution only affects the 
morphology of the SNR at a late time and does not alter significantly the impact on the companion 
star. Thus, all the runs we describe later ignore the binary evolution. The orbital velocity and 
spin, however, are still considered. 

To introduce a SN la explosion, we remove the WD particle cloud and replace it with a high- 
density and high-temperatur e gas. The phys i cal qu antities used to describe the SN la explosion are 



taken from the W7 model bv lNomoto et al.l (jl984l ). It has a mass Msn = 1.378Mq, total explosion 



energy Egn = 1.233 x 10^^ erg, and average ejecta speed Vsn = 8.527 x 10^ km s~^. The kinetic 
energy of the WD from the orbital motion also is added to the SN. We rely on nickel (^^Ni) as a 
tracer for the SN la material. 

Although we have AMR, the WD is still much smaller than the minimum zone spacing. Thus, 
a sub-grid model is required to approximate the SN la explosion. To minimize grid effects, we use 
a small spherical region with a radius equal to fifteen times the smallest zone spacing to represent 
the SN material. Even if we put in the same amount of explosion mass and energy, the results 
are actually sensitive to the initial distribution of the SN ejecta within this small region. The 
detailed setup of the sub-grid SN model and the corresponding behavior will be discussed in the 
next section. 



3. RESULTS 

In this section, the evolution of the companion stars during the SN la explosion for the MS 
star, RG, and He-star binary companions (all the simulation cases are summarized in Table [2|) is 
qualitatively described. Convergence tests of the two-dimensional and three-dimensional simula- 
tions were performed to determine the sensitivity of the numerical spatial resolutions. A detailed 



description of the SN la explosion subgrid model setup is examined as well. 



3.1. A qualitative description of the evolution after the SN la explosion 



Figure [2] demonstrates a typical evolution of the gas density in the orbital plane for the impact 
of a SN la in the MS- WD scenario. In this case (the case MS-Nr in Table [2j), the simulation is three- 
dimensional, but the orbital motion is ignored. The initial binary separation is 1.51 x 10^^ cm, 
~ 3i2*, and 6/8 levels of refinement are used. At 88 sec, the SN ejecta reaches the MS star 
companion and forms a bow shock at the leading surface of the MS star (label A in Figure [2]). The 
bow shock extends further, creating an opening angle with ~ 40° with respect to the x— axis. As 
the shock propagates further, the ejecta sweep around the MS star. At about 200 seconds, the 
ejecta fully surround the MS star and self-interact at the back side of the MS star (label B Figure [2] 
at 203 sec). Subsequently, a cone-shape tail shock forms (label B at 431 sec). The shock waves 
also penetrate the MS star (label C in Figure [2]), reaching the center at about 431 seconds. When 
the shock passes through the MS star, a reverse shock begins to reflect and refill the central SN 
region (label D in Figure [2]). The reverse shock mainly contains companion material and forms a 
solid angle spanning about 40° — 50° with respect to the x— axis. After the impact, the MS star 
pulsates, creating shocks around the surface of the MS star (label E in Figure [2]). By the end of 
the simulation, O.I8M0 of mass has been lost (gravitationally unbound) due to the impact of the 
SN la. The MS star also receives linear momentum and hence a kick of 112 km/sec from the ejecta. 



Qualitatively, the simu lation resembles the results in [Marietta et al.l (J2000l ) and lPakmor et al. 



(J2008l ) . The case rp3_20a in lPakinor et al.l (12008) rese mbles ours, except that the orbital separation 



is smaller (a = 2.68 x 10^^ cm in lPakmor et al.l ( 20081 )). The clearer shock structures reveal typical 



behavior when comparing grid-based hydro and SPH simulatio ns. The most significa nt difference is 
the absence of a reverse shock in the SPH simulations. Although lPakmor et al.l (J2008l ) demonstrated 
that their SPH code could produce the same amount of unbound mass as a grid-based hydro code, 
we believe that the high-speed reverse shock actually plays a more important role for the dynamics 
of the remnant star as it could take away some momentum from the SN ejecta. 



3.2. Effects due to the subgrid SN la setup 



In SNe la, hydrogen has not been detected. However, hydrogen-rich material from the accretion 
disk or from mass stripping during the SN impact is possible and could be detectable when the SN 
ejecta become transparent. Thus, one of the most significant quantities we wish to determine is the 
amount of mass lost from the hydrogen-rich companion star after the SN la explosion. The unbound 
mass could be found by calculating the difference between the initial mass of the companion star 
and the integrated total bound mass from the companion at each timestep. The total bound mass 
is the sum of the companion cell masses in all zones for which the total energy is negative. 



10 




Fig. 2. — Gas density distribution in the orbital plane for a three-dimensional SN la simulation with 
a MS star binary companion (case MS-Nr in Table[2]). The initial orbital separation is 1.51 x 10^^ cm 
(RLOF), but orbital motion is ignored. Each frame shows a portion of the domain spanning 20i?* 
at the labeled simulation time. The color scale indicates the logarithm of the gas density in g cm~'^ 



11 



0.25 




2000 



4000 



6000 



8000 



10000 



0.5 



„ 0.4 



0.3- 



o 

c 
Z3 



0.2 



0.1 



0.0 



MS-WD 





2000 



4000 6000 

Time (sec) 



8000 



10000 



Fig. 3. — Unbound mass vs. simulation time for different SN explosion setups (UPPER: HCV; 
LOWER: MS-WD). CV denotes the case with constant density and radial velocity distribution; LV 
denotes the case in which the velocity is linearly increased in radius, but the density distribution is 
uniform; W7 denotes the fitted model described in section [2^ The simulations are two-dimensional, 
and the binary separation is 2.03 x 10^^ cm (— "^ ° "> ^^^ ^^'^ un-^r^ „„„^„^;„ „^a i ki -.^ mil 
(~ 2.74R^) for the MS-WD scenario. 



3i?*) for the HCV scenario and 1.51 x 10 cm 
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Because of the limitation of spatial resolution, the explosion cannot be perfectly modeled. It is 
found that the amount of unbound mass is sensitive to the sub-grid model representing the SN la. 
Therefore, a realistic description of the SN la explosion is required in addition to including the 
correct amount of energy. 



Dwarkadas &: Chevalier! (119981 ): iDwarkadad (|200[l ) examined the interaction of SN la ejecta 
with the surrounding interstellar medium, suggesting a density profile characterized by an expo- 
nential for setting up the SN la explosion. However, our simulations are restricted to a much smaller 
scale (10^'' cm to 10^^ cm rather than 10^^ cm); the binary companion will break the self-similarity 
of the SN explosion. Thus, the exponential decay may not be appropriate here. We experimented 
with several different initial distributions of SN la in order to study their effect on the post-impact 
companion star. 

The most simple and naive model is to set up a uniform density and constant radial velocity 
distribution in the SN region (denoted by CV SN). However, in reality, the beginning of the ex- 
plosion looks more like a Sedov explosion. Thus, a second model is set up with the radial velocity 
distribution linearly increasing with radius; the gas density is assumed to be uniformly distributed 
(denoted by LV SN). We calculated the amount of unbound mass after the SN la explosion for 
these models and found that the different SN la models will result in a difference of unbound mass 
of about ~ 0.05 - O.SM© (see Figured]). 

Thus, a more realistic sub-grid model is necessary to model the SN la explosion . The early 

time evolution and distribution of the W7 mode l is well descri b ed in iNomoto et al.l (jl984l ) . We 

approximate the density profile in Figure 4b of iNomoto et al.l (j 19841 ) by following a power-law 

distribution with a fixed slope in the M — r relation when the explosion expands to the size of our 

sub-grid model. Re- For a given total mass (Mwd) and an explosion size {Re), the density profile 

can be calculated by 

Mo 



p{r) 



+ CR: 



(2) 



where Mp = 2.4998 x 10^^ g is a constant only related to the slope of the M — r relation in 



Nomoto et al.l (| 19841 ). The other constant, 

c = 



J 47rA/Q 



1 



(3) 



depends on A/q and Mwd- The velocity distribution also is assumed to increase linearly in ra- 
dius, v{r) = fvsn(r /Rr), which is a reasonable assumption and comparable with Figure 4c of 



Nomoto et al.l (| 19841 ) . If we know the density distribution and the total kinetic energy, the coeffi- 
cient of the velocity profile, /, can be calcu lated by integr a ting t he kinetic energy density. Based 
on the temperature profile in Figure 4a of INomoto et al.l ()1984l ). we assume the temperature is 
uniformly distributed. We use the Helmholtz EOS solver to iterate the temperature until the total 
explosion energy is equal to 1.233 x 10^^ erg (denote by W7 SN). 



Figure [3] illustrates two-dimensional simulations of the unbound mass versus simulation time 
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with different SN la models and different companion models (HCV and MS- WD). For the HCV 
companion model, the final unbound mass is O.llS Mp) for the W7 SN mo del, 0.175Mq for the 
CV SN model, and 0.22 5Mf:, for the LV SN model. [Marietta et aD (|20od l obtained O.ISM© by 
using their HCV model; iPakmor et al.l (|2008l ) using the same model but with an SPH simulation, 
obtained O.134M0. Their final unbound mass is between our CV SN model and W7 SN model. 

For the other companion model, the MS- WD scenario, we obtained O.l JdMp) for the W7 SN 
model, O.33OM0 for the CV S N model, and 0A2 7 Mcr^ f or the LV SN model. IPakmor et aD j^OoJ) 
used t he same compa i iion f rom llvanova &: TaamI (|2004l ) ; due to the difference of constructing this 
model, IPakmor et al.l (|2008l ) has a slightly larger radius and different material compositions. From 
their fitted power-law relation, the final unbound mass at our initial binary separation (a = 2.74/?*) 
is O.I9M0, which is consistent with our results. We therefore use the W7 SN model for all other 
cases in this paper. 



3.3. Convergence Test 

In order to check the sensitivity of unbound mass to spatial resolution, we performed con- 
vergence tests on two-dimensional and three-dimensional simulations using the MS- WD com- 
panion model. We performed several simulations with the same initial binary separation (a = 
1.51 X 10^*^ cm), and the same W7 SN model for the SN la explosion, but with several different 
AMR levels, computing the amount of unbound mass as a function of time for each simulation. In 
Figure m we see that the unbound mass decreases when the zone spacing decreases; the decrement 
changes very little from 6/8 to 7/9 levels. The difference between 6/8 (corresponding to 1024^ in 
3D and 1024 x 2048 in 2D) and 7/9 (corresponding to 2048^ in 3D and 2048 x 4096 in 2D) levels 
is within 3%. The final unbound mass of this two-dimensional run for the MS- WD scenario is 
Mub = O.175M0. 

For the convergence test of three-dimensional simulations, we use the same companion model 
(MS- WD), but the spin and orbital motion are included. In this stage, we are only interested in 
the sensitivity to numerical resolution during the SN la explosion. Thus, we ignored the effects 
of accretion of Roche-filling material, i.e. the WD exploded immediately when it was added to 
the simulation. The result is shown in Figure [H The convergence behavior is the same as for 
two-dimensional simulations. However, three-dimensional runs require a higher AMR level for the 
same spatial resolution because the cylindrical mesh has only half of the zone spacing, due to 
symmetry, of the Cartesian mesh for the same number of AMR levels. The final unbound mass 
Mub = 0-195Mq for 7/9 levels. The difference between 6/8 and 7/9 is about 8%, which is larger 
than the difference in the two-dimensional case. However, a full set of runs at 7/9 levels proved to 
be too expensive for us to undertake. Thus, we chose 6/8 levels for all other production runs. 
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Fig. 4. — Unbound mass from the companion star vs. simulation time using different AMR re- 
finement levels for the MS- WD scenario. The initial binary separation is 1.51 x 10^^ cm and the 
calculation is two-dimensional. 
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Fig. 5. — Similar to FigurcHJ but with three-dimensional simulations, including orbital motion and 
spin. 
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3.4. Evolution of the Companion Star after the SN la Explosion 

We have performed several simulations with the MS- WD, RG-WD, and He-WDc companions 
(see Tabled]). For each scenario, we ran cases with and without the orbital motion to distinguish the 
effects of asymmetry. We also compared simulations between three dimensions and two dimensions. 
In addition, a parameter survey of changing initial orbital separations was also carried out. 

Table [5] summarizes all our numerical simulations. The cases beginning with "MS," "RG," 
and "He" correspond to the MS- WD, RG-WD, and He-WDc scenarios in Table [H respectively. 
The letter "r" represents the cases with initial orbital separation equal to the separation for RLOF 
(see equation [1]) ; "4" and "5" are cases with initial binary separations equal to 4i?* and 5-R* . The 
corresponding orbital speeds are calculated for initially circular orbits using Kepler's third law. 
The letter "N" represents cases that ignored the orbital motion and spin. The cases with "2D" are 
two-dimensional simulations; others are simulated in three dimensions. 



3. 4-1- Main-sequence binary companion 

In Section 13.11 we described the case MS-Nr (see Table [5]) for a three-dimensional simulation 
with the MS- WD companion but without orbital motion or spin. Figure [6] shows a similar simulation 
for the case MS-r, which includes the orbital motion and spin. The simulation resembles the case 
MS-Nr during the first few hundred seconds, except that the tail shock is shifted and curved because 
the ejecta speed is much higher than the orbital speed. However, after a thousand seconds, the 
asymmetry becomes more obvious; most of the same features, however, can still be seen. This 
asymmetry will become more important for the evolution of the SNR. Including orbital motion and 
spin leads to a greater unbound mass by about 16%, but the kick velocity remains the same. The 
additional unbound mass may be due to the spin; the rotational energy makes fluid elements easier 
to unbind. However, the momentum lost in this additional mass may not change the center-of-mass 
velocity. A more detailed description and explanation is given in Section [4.11 



3.4-2. Red- giant binary companion 

The evolution of gas density for the RG-WD companion model is shown in Figure [7] (RG- 
Nr, without orbital motion) and Figure [8] (RG-r, with orbital motion). The asymmetric effect is 
insignificant for the RG companion, because the orbital speed is much lower than the ejecta speed. 

Unlike the MS star, almost all the envelope (> 95%) of the RG is removed during the initial 
impact. Although the ejecta speed is smaller than for the MS- WD companion when it reaches the 
RG, the greater extent and lower binding energy of the RG makes it difficult to survive the SN la 
explosion. As a result, most of the momentum from the SN ejecta is transferred to the unbound 
mass, creating only a tiny kick on the RG core (less than 30 km/sec). The bow shock has an 
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Fig. 6. — Similar to Figure O but for case MS-r; each frame reveals a portion of the domain 
spanning 30i?*. 
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Fig. 7. — Similar to Figured! but for case RG-Nr. 
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Fig. 8. — Similar to Figured! but for case RG-r. 
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opening angle of about 40° — 45° with respect to the x— axis. For the RG-r case, a spiral pattern 
appears due to the combination of fallback and spin about 5 days after the SN la explosion. The 
reverse shock also is weaker than in the MS- WD scenario. We can only report a lower limit of 
the final unbound mass: the core of the RG is represented by a particle cloud and the mass of the 
particle cloud is not lost in our simulation. 



3.4-3. Helium star binary companion 

The He-WDc companion has the most dramatic impact and is the most compact companion 
in comparison to the MS- WD and RG-WD cases. Figure [9] and Figure [10] show the evolution 
of gas density in the orbital plane for the He-Nr and He-r cases. The opening angle of the bow 
shock is about 30° — 35°. Because the orbital speed is higher than in the MS- WD and RG-WD 
scenarios, asymmetric effects due to the orbital motion appear earlier (within the first few hundreds 
of seconds). 

After the initial impact, the compressed helium star oscillates. The oscillation frequency is 
higher and the flow is more turbulent in comparison to the MS companion. After several hundreds 
of seconds, the strong velocity shear behind the helium star leads to Kelvin-Helmholtz instability. 
Later, the oscillation of the helium star generates shocks which interact with Rayleigh- Taylor 
instabilities owing to the gravity from the helium star. Therefore, at the end of the simulation, 
a turbulent remnant environment is obtained. Figure [11] shows the velocity vector field in a gas 
density slice at 2231 sec in the orbital plane for the case He-r. There is a global kick velocity in the 
lower-left direction, caused by the orbital speed and a random turbulent velocity. 



4. Discussion 

4.1. Parameter Survey 

Although most SDS channels assume the binary system is in RLOF, the actual evolutionary 
stage and explosion conditions may differ from the model we used in our simulations. We therefore 
have performed a parameter survey to explore the dependence of the numerical results on varying 
binary separation. Significant physical quantities for each run are summarized in Table [21 

Figure [12] and Figure [13] show the final unbound mass and companion kick velocity as functions 
of the binary separation in units of companion radius. Note that the unbound masses for the RG- 
WD scenario are only lower limits, because the unbound mass is still increasing at the end of these 
simulations. In addition, the core of the RG is represented by a rigid particle cloud for which 
the mass is assumed to be unchanged during the simulation. Therefore, we effectively ignored 
interactions such as mass stripping and ablation on the core of the RG. The companion kick 
velocity can be derived from the difference between the companion center-of-mass velocity at the 
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Fig. 9. — Similar to Figure [U but for case He-Nr. 
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Fig. 10. — Similar to Figure [U but for case He-r. 
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Fig. 11. — Velocity vector field around the helium star companion at 2231 sec after the SN la 
explosion in the simulation He-r. The velocity vectors are plotted every 10 zones. The white 
arrows reveal the scale of the domain and the scale of the velocity vectors. The color displays the 
logarithm of the gas density in g cm~^ in the orbital plane. 



24 



10" 



10 



•I , 



I 



■•= 10 * 



ID'- 



='«^— ^ 


*^ -1 , ^ \ 








"^^^^K — ^ — ' 


■^^^ 




\ 


\ 


^^^^^^*— ^ 


^-J;^^^^^ 


\ 


'%. 




Lower iinrit 


V 


k ^v] 


^ 




^ 


"^L 1 


^A 






V \ 


^5v 


■ 




> "X 


■ v\^ 






\ -\ 


\v 






s"" Ai 


^^ 






^■'Nu 


\^t 




r 


v.^v 


^^v 


- 




■*■■. X 


^t 








*\ * 






^ \ 








\.x 


^A 






^\ 


Sk 






V X 








'^ %. 


*^x 








\ '\ 








"^ ' vH. """IcL. 


* 






k C\ 








'^L 


A. 






' '^\ 


'^^ 






^L 


^^\,- 








V -^ 






MS 




*^ 


'A. 


■ 


He 




''v^^ 


, 'X 




RG 




\^ 


\, *\ 


* * 


Marietta (00) 






\\ ' 


A A 


Pakmor f08.) 






* ^V 








^v 



10" 



IQ^ 



Sepia rsCFon iit.l 



Fig. 12. — Final unbound stellar mass vs. binary separation for different companion models. Each 
color represents data for a type of companion model: black represents simulations for the MS- WD 
scenario; blue, the He-WDc scenario; and red, the RG-WD scenario. The separation is expressed 
in units of the companion radius. Lines show fitted power-law relations based on the numerical 
simulation results (the red line shows the lower limit of the RG-WD scenario). The "+" symbols and 
dashed lines indicate 3D data without orbital motion; "o" symbols and solid lines present 3D data 
with orbital motion; and "x" symbols and dashed dot line s demo nstrate 2D data without orbital 
motion. Star symbols show the data from Marietta et al.l (J2000) for the HCV mode l. Triangles 
indicate the power-law relation for a MS binary companion from lPakmor et al.l ([200a). 
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Fig. 13. — Similar to Figure [T2l but for kick velocity vs. binary separation in different companion 
models. 
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end of a simulation and the velocity at the time when the ejecta reached the companion star in 
that simulation. 

In general, the final unbound mass and kick velocity decrease when increasing the binary 
separation. We find that the final unbound mass, 5Muh, can be described by the relation 

6AU = C^^a'^-^MQ, (4) 

where a is the orbital separation, ?7T,ub is the power-law index, and the constant Cub depends on the 
companion model (see Table [3]). All the companion models we used could be described with this 
power-law relation; the power-law indices, however, are different for different companion models. 

Although the initial binary separation is the smallest for the He- WD scenario, suggesting the 
impact of the SN ejecta for this scenario should be the largest, the He star is more compact and 
has a higher binding energy than the MS star or RG. Therefore the amount of unbound mass is 
smallest in the He- WD scenario. In the RG-WD scenario, the gravitational force acting on the core 
region is only approximate since the softening length in the particle cloud representing the core is 
comparable to the size of the remaining envelope. However, most of the envelope is lost with at 
least 75% — 96% ejected depending on the binary separation. 

For the MS- WD scenario, the difference between two-dimensional and three-dimensional runs 
varies from 3%, for the smallest binary separations, to 10%, for the largest binary separations. The 
error in the unbound mass is usually higher for the larger binary separations because the unbound 
mass decreases by a power-law relation with a high power-law index. The results between 2 and 3 
dimensions in the He- WD scenario are even larger than in the MS- WD scenario. The uncertainty 
arises not only from the low percentage of unbound mass for the He star, but also from the high 



amount of turbulence seen in these models. iPan et al.l (|201Cl ) demonstrated that the effect of the 



turbulence precludes convergence of the unbound mass no matter the spatial resolution and leads 
to a 14% error for the He- WD scenario in the two-dimensional simulations. 

By including the orbital motion and spin, the MS- WD scenario loses 16% more mass for the 
MS-r run and 4% for the MS-5 run. The additional mass lost may not be accurate for runs with 
large separations: the mass loss decreases dramatically at larger separations. Other uncertainties 
associated with the spatial resolution (or the gravitational potential in the binary system) may 
be comparable to the percentage of mass difference in runs characterized by large separations. 
In general the orbital motion should not change the binding energy significantly if the binary 
separation is larger than the Roche separation. In addition, the orbital speed is much smaller 
than the SN ejecta speed. Thus, orbital motion should not have much effect on the impact of the 
explosion, instead simply providing a constant velocity for the companion star. However, the spin 
of the companion may have a stronger infiuence. Adding rotational energy to the companion star 
should make its envelope's binding energy lower. This addition may explain the extra mass lost in 
our simulations that include spin. Figure O shows the distribution of the ratio of rotational energy 
density and binding energy density in the x — z plane. The rotational energy density in the back 
and front region of the companion star is ~ 3% of the binding energy density at i? ~ 0.8i?*. When 
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Fig. 14. — The distribution of the ratio of the rotational energy density and binding energy density 
in the x — z plane for the MS-r run. The region is rescaled to the units of the companion star 
(i?* = 5.51 X 10^° cm). 
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Fig. 15. — Initial position of passive particles in the x — z plane. Particles are plotted within a 
slice of thickness equal to three zones. The axes are scaled to the radius of the companion star 
(-Rms = 5.51 X 10^ km). Colors represent the gravitational binding status at t = 144 sec (upper) 
and t = 1060 sec (lower) for the case MS-r (see Table [2]): blue for bound; red for ablated; green for 
stripped; and black for particles already gone from the simulation box at that time. The black grids 
between particles reveal some numerical artifacts in particle initialization, but this artifact does 
not alter the dynamics of the fluid or the tracer particles. A detailed explanation and calibration 
are given in Section [4.21 
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compared with the unbound region in Figure [T5l there is evidence of some overlap of the unbound 
mass with the high rotational-energy density region. Therefore, the spin of the companion should 
broaden the unbound mass region. 



The simulation data from [Marietta et al.l (12000 ) and lPakmor et al.l (120081) a r e also plotted and 
compared in Figure [121 Although the MS companion model in [Marietta et al.l ( 2000l) is different 



from that used here, the general behavior of their results is consistent with ours. iPakmor et al 
(J2008l ) have a MS star model that resembles ours but is slightly larger. Th e power-law relation for 
the final unbound mass versus separation provided bv lPakmor et al.l (|2008l ) also is consistent with 
our MS- WD scenario. 



Leonard (J2007l ) studied the spectra of two SNe la (SN 2005am and SN 2005cf) in order to search 
for the Ha emission that may emerge from the unbound hydrogen-ri ch companion rnass. He relied 
on the data from the two SNe and coupled the modeling re s ults o f iMattila et al.l (J2005l ) and the 
final unbound mass in [Marietta et al.l (J200d ) and lMeng et al.l (J2007l ). providing an estimated upper 
limit of 0.02 Mq for the amount of solar abundance material that may have remained undetected. 
In contrast, the mass loss for our cases in RLOF are all greater than these limits, except for the 
He-WDc scenario. However, the limits should depend on the actual evolutionary stage of the 
companion candidate. If the upper limit is correct, SN 2005am and SN 2005cf may result from a 
more compact binary companion or from the DDS. 

Similarly, the kick velocity versus binary separation can also be described by a power-law 
relation, 

•^^kick = Ckicka'"'^'*M0, (5) 

where mkick is the power-law index, and c^ick is a constant depending only on the companion model 
(see Figure [T3|) . The power-law indices and constants for the final unbound mass and kick velocity 
are summarized in Table [31 For the RG-WD scenario, because the simulation domain is of a very 
large scale, the error in kick velocity is comparable to the value we deduced. Therefore, we can 
only provide an upper limit on the kick velocity of about 40 km/sec. 

If we include the orbital motion and spin, they only produce about a 2% difference in the 
kick velocity, which suggests that the orbital motion and spin are unimportant to the momentum 
transfer between ejecta and companion. 

Although our final unbound mas s is consistent with Pakmor et al.l (120081 ) . the kick velocity we 
determined is lower than the kick in IPakmor et al.l ( 20081 ) . This result may be due to the lack of 
a reverse shock in the SPH simulations: the reverse shock carries some momentum away from the 
SN ejecta, consequently making the companion kick smaller. 
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Table 1: The non-degenerate companion models. 



Model 



(Mo) 



(km) 



[g/cvo?) 



pa 

(dyne/cm^ 



Ti c 
c 

(K) 



(sec) 



HCV 


1.027 


6.76 X 10^ 


112.3 


1.90 X 10^^ 


1.48 X 10^ 


1.54 X 10^ 


MS- WD 


1.17 


5.51 X 10^ 


68.9 


1.33 X 10^^ 


1.44 X 10^ 


1.06 X 10^ 


RG-WD 


0.98 


2.15 X 10^ 


3.17 X 10-3 t 


5.50 X 10"t 


1.88 X W^^ 


2.82 X 10^ 


He-WDc 


1.007 


1.36 X 10^ 


6.10 X 10^ 


5.12 X 10^9 


1.28 X 10^ 


1.40 X 10^ 

















"Stellar mass, M 

Stellar radius, i?, 
'^Central density,pc 
''Central pressure, Pc 
'^Central Temperature, Tc 

■'^The dynamical time scale, tdyn = l/2(Gp)^^'^ 
^For the RG-WD model, the central density, pressure and temperature are the central values of the envelope mapped 

in FLASH, not the core values in MESA. 




-♦*■ 



Fig. 16. — Similar to Figure \T5\ but for case He-r at i = 265 second. 
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Table 2: Simulations 





a 


a 


JM^ftb 


SVkick "" 


Vorb '^ 


Mstr/Mabl ^ 


Case 


(cm) 


(Mo) 


(km/sec) 


(km/sec) 


(Mo) 


MS-r 


1.51 X 


10" 


2.08 X 


10-1 


114 


256 


0.65 


MS-4 


2.20 X 


10" 


5.00 X 


10-2 


58.8 


212 


0.43 


MS-5 


2.75 X 


10" 


2.41 X 


10-2 


39.6 


190 


0.73 


MS-Nr 


1.51 X 


10" 


1.80 X 


10-1 


112 





0.49 


MS-N4 


2.20 X 


10" 


4.35 X 


10-2 


57 





0.39 


MS-N5 


2.75 X 


10" 


2.31 X 


10-2 


41 





0.58 


MS-2D-Nr 


1.51 X 


10" 


1.84 X 


10-1 


113 


— 


— 


MS-2D-N4 


2.20 X 


10" 


4.21 X 


10-2 


59.5 


— 


— 


MS-2D-N5 


2.75 X 


10" 


2.08 X 


10-2 


36 


— 


— 


RG-r 


6.13 X 


10^2 


6.37 X 


10-1 


12 


41.8 


0.01 


RG-4 


8.58 X 


1012 


5.93 X 


10-1 


25 


35.3 


0.06 


RG-5 


1.07 X 


1013 


5.07 X 


10-1 


19 


31.6 


0.19 


RG-Nr 


6.13 X 


1012 


6.37 X 


10-1 


14 





0.01 


RG-N4 


8.58 X 


1012 


5.97 X 


10-1 


23 





0.046 


RG-N5 


1.07 X 


1013 


4.99 X 


10-1 


4.6 





0.16 


He-r 


3.86 X 


lOio 


2.3 X 


10-2 


65.6 


522 


0.82 


He-4 


5.43 X 


lOio 


8.1 X 


10-3 


28.9 


440 


0.74 


He-5 


6.79 X 


lOio 


4.11 X 


10-3 


14.7 


394 


0.74 


He-Nr 


3.86 X 


lOio 


1.99 X 


10-2 


77 





0.48 


He-N4 


5.43 X 


lOio 


7.3 X 


10-3 


29 





0.53 


He-N5 


6.79 X 


lOio 


3.82 X 


10-3 


15.5 





0.56 


He-2D-Nr 


3.86 X 


lOio 


1.73 X 


10-2 


63 


— 


— 


He-2D-N4 


5.43 X 


lOio 


5.3 X 


10-3 


12.7 


— 


— 


He-2D-N5 


6.79 X 


lOio 


2.6 X 


10-3 


6.3 


— 


— 



















"The initial orbital separation, a 

The final unbound mass, SAIuh, described in equation ID 
''The kick velocity, Svkick, described in equation [5l 
■^The initial orbital speed 
''The stripped mass to ablated mass ratio 
^Note: The cases beginning with "MS," "RG," and "He" correspond to the MS- WD, RG-WD, and He-WDc scenarios 

in Table [T] respectively. The detailed description for the notations of simulation cases is described in Section 13.41 

Orbital motion and passive particles are not included in the two-dimensional simulations, so for these cases the orbital 

velocity and stripped to ablated ratio are not calculated. 



32 



ja 



•A 






-** 



•^m 




* • _ » .■»* 



^*5r 






oa 



T.JI 



Fig. 17. — Similar to Figure [T5l but for case RG-r at t = 16.0 hours. 
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4.2. Stripped and Ablated Mass 



The unbound mass could be either due to ablation (heating) or to stripping (mo mentum trans- 



fer), depending on the energy of the SN ejecta and the binding energy of the companion !) Wheeler et al 



19751 ). Quantitatively determining the ablated and stripped mass requires a method to trace the 
fluid elements in Lagrangian coordinates. This method can be implemented in FLASH using pas- 
sive particles. We create an additional particle property for each passive particle and then label 
its gravitational-binding status as bound, ablated, or stripped. 10^ passive particles are randomly 
distributed, with particle density proportional to the gas density. All passive particles are initial- 
ized as bound and then the gravitational-binding status in each timestep is recalculated. Particles 
in zones with negative total energy (-Etot = -E'grav + -E'kin + -^'thermal) are considered to be bound; 
otherwise, we mark passive particles in zones as ablated (stripped) if the thermal energy is larger 
(smaller) than the kinetic energy. Ablated and stripped particles possibly become bound again in 
the subsequent timesteps, but they are not allowed to switch from one state to the other. 

Figure [15] illustrates the initial distribution of passive particles colored by their gravitational- 
binding status in the x — z plane. At 144 seconds (upper panel), particles on the front side of the 
companion are all ablated as a result of the initial impact. The MS star companion is compressed 
and ablated to the initial radius of ~ 0.4i?*. This result corresponds to the first peak of the unbound 
mass (~ 0.2Mq) in comparison to the simulation time in Figure[5j Next, the reverse shock, traveling 
in the positive x— direction, thins the unbound layer to ^ 0.5i?* and decreases the unbound mass. 
Following this decrease, the companion is surrounded by the ejecta while mass around the surface 
and in the tail region becomes unbound again (lower panel at 1060 sec), eventually increasing to 
O.2IM0 at the end of the simulatio n (Figure [5l). This non -monotonic behavior of changing unbound 



mass in time also has been seen by lPakmor et al.l (|2008l ) in their SPH simulations. They, however. 



do not examine the unbinding process at early stage and their non-monotonic behavior is relatively 
weak. Therefore, they deemed it to be a numerical artifact in their convergence study. However, it 
is difficult to definitively show that the non-monotonic behavior is numerical in origin. 

In contrast to the MS-WD scenario. Figure [16] and Figure [T71 reveal similar distributions for the 
He- WD scenario and the RG-WD scenario. The He- WD scenario's behavior resembles the MS-WD 
scenario, but mass is only ablated in the leading surface region. More stripped mass in the tail 
region results from the high turbulence around the He star. The dramatic impact and turbulence 
make momentum transport more efficient. In the RG-WD scenario, almost all the envelope is 
ablated during the initial impact; only a small amount of mass (~ O.OOTMq) is stripped. Later 
on, some mass around the core region becomes bound again after the RG is compressed. After the 
impact, the bound mass slowly increases because of fallback. 

The amount of ablated and stripped mass can be calculated by counting the ablated and 
stripped particles. Non-interacting passive particles do not affect the accuracy of the subsequent 
tiine evolution because passive particles simply move in Lagrangian coordinates without interact- 
ing with fluids. However, the initialization of passive particles in our setup was found after our 



34 



simulations to have some numerical artifact related to block boundaries. This grid-shaped artifact 
can be seen in Figure [T5l Therefore, a weighting factor for the ablated/stripped mass due to each 
particle is necessary when calculating the total particle mass. 

We calculate the particle weights by evaluating the difference between particle density and gas 
density in each zone for the initial conditions. Thus, the particles in any given zone have the same 
weight. However, in low-density regions, some zones may not contain any particles; this lack may 
lead to another uncertainty, underestimating the total particle mass. Thus, we add this missing 
mass by assuming that the stripped mass to ablated mass ratio in these zones is the same as for 
their neighbors. After this correction, we are able to calculate the stripped mass and ablated mass 
for each timestep. 

The stripped mass to ablated mass ratios for all the runs are summarized in Table El These 
values are calculated from the unbound particles at a time after the initial impact and before the 
unbound particles leave the simulation box. It is found that the unbound mass is mainly due to 
ablation in all the runs. The stripped mass to ablated mass ratio is about 0.4 — 0.7 for the MS- 
WD scenario, 0.5 — 0.8 for the He- WD scenario, and less than 0.2 for the RG-WD scenario. This 
result contradict s our previous results using the mixing of SN ejecta at late times in 2D simulations 
(jPan et al.l l2010l ) (tracer particles are not supported in 2D). The ablation happens at the very 
beginning of the initial impact, but the mixing of ejecta occurs later, based on our analysis using 
passive particles. Thus, the mass treated as stripped in the previous study may already have been 
ablated during the initial impact. In addition, abla tion was previously ignored in most previous 
analytical or semi-analytical work (JMeng et al.l 120071 ) , resulting in an underestimate of the impact 
of SN ejecta on the binary companion. 



4.3. Hole in the Supernova Remnant 

In the SDS, the impact of SN ejecta not only can affect the companion star, but it also can 
affect the shape of the ejecta: a cone-shaped hole shadowed by the companion will break the sym- 
metry of the ejecta. Observationally, spectropolarimetry studies reveal a contin uum polarization of 
about 0.2% - 0.7% for normal SNe la and 0.3% - 0.8% for subluminous SNe la ( Howeh et allboOl 



Wang &: Wheeleiil2008l ). implying that the asymmetry in SN la explosions is small but real. Recent 
theoretical and numerical studies suggest that the intrinsic variation in SN la light curves poten- 
tially results from yiewi r ig asymmetric explo sions from diffe rent angles (Kaseri et al. l2004l : iKasen 
2010l : iMaeda et aUboiol : blondin etaPboilh . Furthermore, ICarcia-Senz et al.l (|2012l ) studied the 
interaction of the hole, SN material, and ambient medium using SPH simulations. They conclude 
that the hole could remain open in the SNR for hundreds of years, suggesting the hole could affect 
its structure and evolution. 

Figure [18] shows the companion-material and SN-material density distribution in the orbital 
plane for the MS-r, RG-r, and He-r cases. The companion material is confined mainly in the bow 
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Fig. 18. — Companion and SN gas density in the orbital plane for different scenarios. The left 
panel displays the density of the companion material (^ff +^ He +^^ C) and the right panel shows 
the density of SN material (^^Ni). Companion model and initial setup from the top row to the 
lowest row are MS-r, RG-r, He-r cases in Table [2j 
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Fig. 19. — Speed distribution of unbound mass and compositions for the MS-r, RG-r, He-r cases 
(Table [2|). Each color represents the speed distribution for a particular run: black represents the 
case MS-r; red, case RG-r; and blue, case He-r. 
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shock; the opening angle of the hole is about 40° — 45° for the MS-r case, 45° — 50° for the RG-r 
case, and 30° — 35° for the He-r case. The MS star and RG companions have a similar opening 
angle, but the high-speed tail shock in the RG companion has the largest opening angle due to 
having the lowest ejecta Mach number. The high orbital speed and strong turbulence in the He-r 
case distort and mix the hole with more SN material than other companion models, creating a 
smaller hole with more fallback in the SNR. 

We may also ask how the unbound hydrogen-rich material is distributed in the SNR and how 
this material can be detected by observations. Based on our simulations, the unbound mass is 
mostly confined in the hole. The velocity distribution of unbound mass is shown in Figure [T9l The 
peak velocity is 550km/sec~^ for the MS- WD, 955km/sec~^ for the He-WDc, and 660km/sec~^ 
for the RG-WD. The RG companion has the largest velocity dispersion, but the peak velocity 
is slightly larger than the peak velocity for the MS star companion. This low outflow velocity 
(< l,000km/s"^) is much smaller than the ejecta speed (~ 8,000km/s~^), suggesting that the 
hydrogen could be hidden in the ejecta in the early-stage explosion, except when looking directly 
from the backside of the companion. However, it is still possible to detect it at late times when the 
expansion slows down and the ejecta becomes transparent. 



4.4. Nickel Contamination 

The companion star could be contaminated by the SN ejecta in its envelope during the initial 
impact or as fallback at late times. If the mixing of SN ejecta occurs mainly in the companion's 
envelope, the contamination provides one possibility for detecting iron abundance from the remnant 
star after the nickel radioactively decays. Since we use nickel as a tracer for the SN la ejecta in our 
simulations, we could estimate the nickel contamination in the remnant star in a SNR by calculating 
the amount of bound nickel after the SN la explosion. Following the same method as that used to 
determine the bound mass, we can calculate the amount of bound nickel by summing the ejecta 
material that cannot escape the gravitational potential. There are essentially two bound regions in 
the simulations: the region around the post-impact companion star and the region at the explosion 
center, separated by the reverse shock. The bound nickel in the SN region likely results from the 
fact that the some of the nickel is bound to the system in the initial setup due to the presence of 
the companion star. We assume that the bound nickel in the companion star more likely reveals 
the nickel contamination. 

Figure [20l illustrates the amount of bound nickel in these two regions for the case MS-r. The 
bound nickel is dominated by the SN region at the beginning and then reaches a peak of ~ 1O~^M0 
at about ~ 80 sec. However, after the interaction between the reverse shock and the SN region, the 
bound nickel in the SN region becomes unbound and eventually disappears after around 2, 000 sec. 
The remaining bound nickel in the MS star companion is ~ 10~^Mq. In general, increasing the 
initial binary separation leads to a greater contamination because a slower ejecta speed makes 
it more difficult for ejecta material to escape; a larger separation, however, leads to a smaller 
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Fig. 20. — Bound nickel around the MS star companion region and the SN region. The plus sign 
indicates the companion region, and triangles indicate the SN region. Total bound nickel is plotted 
using the solid line. 
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companion cross-section, thereby blocking less SN material. Thus, a simple variation for the nickel 
contamination with the binary separation is not found. All the results are consistent to within an 
order of magnitude with a value of ~ W~^Mq. Orbital motion and spin only introduce a minimal 
effect. This amount of nickel contamination is significantly less than the previous estimate by 



Marietta et al.l ( 20001 ) for the HCV scenario (an upper limit of 10 ^Mq) because the interaction of 



the reverse shock. 

Figure [21] shows the nickel contamination for the RG companion, the RG-r case. In this 
example, the bound nickel in the SN region only exists for a very short time at about 10^ sec after 
the explosion. There are two bumps of bound nickel around the RG region. The first bump, at 
about ~ 4, 000 sec, is related to the initial impact of the SN la explosion. Bound nickel reaches the 
first peak when the shock compresses the leading surface of the RG, and then it decreases when 
the reverse shock passes through it. The second bump occurs later at about ~ 10^ sec in a region 
at the back side of the RG while the compressed RG relaxes. The nickel contamination for the case 

0- 



RG-r is about 10"*^ ~ 5 x 10"^ M 



Similarly, Figure [22] shows the nickel contamination for the helium star companion (case He-r). 
The nickel contamination in the He star region stays around ^ 1O~^M0 for the whole simulation, 
but in the SN region it decreases from 10~^ ~ 1O~^M0 to ~ 3 x 10~^. Unlike the MS star cases, the 
bound nickel in the SN region does not disappear at the end of the simulation while interacting with 
the reverse shock. It is more likely that the bound nickel in the SN region will merge eventually 
with the He star region at late times, resulting in a total contamination of ~ 4 x 10~^Mq. Figure [23] 
shows the bound nickel density in the orbital plane at 2, 231 sec for these two regions; at that time, 
the two regions are already mixed with each other. 

If the nickel contamination is restricted to the envelope of the companion star, we ca n estimate 



the n ickel to hydrogen plus helium ratio by using the bound nickel and the envelope mass (jPan et al. 



2010l ). We define the envelope radius using the extreme in second derivative of the companion 
gas density with respect to the radius in the initial conditions. Thus, the envelope mass is the 
mass outside this envelope radius. If we assume the unbound mass is entirely from the envelope 
(Figure [T5]) . then the final envelope mass is the difference between the initial envelope mass and 
the final unbound mass in Table [2] However, we note that the envelope mass (~ O.3M0 it O.IMq) 
we estimate for the MS star companion is close to the final unbound mass (~ 0.2Mq); thus for 
this case we only can provide an order of magnitude estimate. The estimated upper limit for the 
ratio of nickel to hydrogen plus helium is about 10~^ for the MS- WD, 10~^ for th e He-W Dc, and 
2 X 10~^ — 10~^ for the RG-WD scenario. The solar ratio from [Anders fc Grevessd (|l989l ) is about 



5 X 10 ^, which is similar to or slightly smaller than the results in our simulation values, suggesting 
a possible probe to identify the progenitor candidate in a SN la SNR. 
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Fig. 21. — Similar to Figure [20l but for RG binary companion. Note that there is only one triangle 
data point for SN region because the duration of this region is short. We have only one checkpoint 
file during that time. 
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Fig. 22. — Similar to Figure [20} but for helium star binary companion. 
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Table 3: Power-law indices for the final unbound mass and kick velocity. 



Model 


m^h "" 


Cub" 


mMck" 


^kick 


MS 


-3.61 


0.894 


-1.76 


2.83 


MS-N 


-3.46 


0.757 


-1.69 


2.78 


MS-2D-N 


-3.66 


0.859 


-1.88 


2.89 


RG 


-0.391 


-0.0103 


- 


- 


RG-N 


-0.416 


0.00269 


— 


— 


He 


-3.04 


-0.260 


-2.62 


3.02 


He-N 


-2.92 


-0.381 


-2.83 


3.17 


He-2D-N 


-3.36 


-0.245 


-4.6 


3.9 



"The power-law index, rriuh, and constant, Cub, described in equation]!] 
''The power-law index, rrikick, and constant, Ckick, described in equation [S] 
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Fig. 23. — Bound nickel density distribution in the orbital plane around the helium star companion 
(left) and the SN region (right). The black line describes the boundary between the companion 
material region and the SN material region. The color scale indicates the logarithm of the nickel 
density in g cm~^. 
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Fig. 24. — Angular momentum versus simulation time for the MS-r case. The angular momentum 
is calculated in the companion center of mass frame and only considers the bound companion gas. 
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4.5. The Remnant Companion Star 

At the end of our simulations, all the companion stars have survived the impact of the SN la 
explosion. The post-impact velocity of the remnant star will be the original orbital velocity plus 
the kick velocity. In general, this velocity should be larger than the velocities of background stars. 
In our simulations, the kick velocity contributes up to 45% of the final net velocity for the MS-r 
case, 13% for the He-r case and < 50% for the RG-r case. Therefore if the binary systems in the 
SDS are in RLOF, then the final velocity including kick and orbital speed will provide information 
in the search for the progenitor star in the SNR. 



Ruiz-Lapuente et al.l (J2004l ) observed stars near the remnant center in Tycho's SN. They found 
a star named Tycho G which has a higher radial velocity (108 it 6km s~^, in the local standard 
of rest) than other background stars at the corresponding distance (20 — 40 km s~^ at 2 — 4 kpc), 
suggesting Tycho G could be the non-degenerate remnant star resulting from the SDS. Tycho G is a 
slightly evolved star with a stellar type G0-G2, mass M ^ IMc^, r adius i? ~ 1 — 3i?0, and effective 
temperature Tefj = 5, 900±100 K (JGonzalez Hernandez et al.ll2009l ). Thus, if it is indeed a remnant, 
the progenitor star was likely a MS star. Although our simulation is not set up for Tycho's SN, the 
MS-r case resembles it. The final net velocity of our MS-r case is 277km s~^, which covers the value 
of Tycho G. If the remnant star in Tycho's SN has not changed significantly since the SN explosion, 
the property of our post-impact companion may provide further information for identifying the 
progenitor star of Tycho's SN. The post-impact central density, temperature, and pressure in our 
simulation for the MS-WD n iodel are 24 g cm~^, 7 .8 x 10^ K, and 2.7 x 10^^ dyne cm~^, which are 
consistent with the values in [Marietta et al.l (|2000l ). 



As iMeng &: Yang) (120111 ) pointed out, the rotational velocity of the remnant star could be an 
important diagnostic for the Tycho G star as the companion candidate in Tycho's SN with the 
SDS. In their model, the mass, space velocity, r adius, luminosity, and effect ive t emperature of the 



remna nt star are consistent with observations bv iRuiz-Lapuente et al.l (120041) and iKerzendorf et al 



(|2009l ). but not the rotational velocity. The rotational velocit y observed bv IKerzendorf et al.l (J2009|) 
is only 7. 5ifc2 km s~^, but the predicted rotational velocity in lMeng &: Yana (J201ll ) is ~ 100 km s~^. 
However, iMeng &: Yang (120111 ) ignore the impact of the SN ejecta on the companion star and assume 
the remnant star did not change at the onset of the SN la explosion. Our simulation shows that 
about 18% of the mass is lost in the MS-r case, suggesting that the angular momentum of the 
remnant star also should decrease after the impact. 

Figure [Ml shows the change of angular momentum versus simulation time in the companion 
center of mass frame. After the impact, the MS companion loses 48% of its initial angular momen- 
tum (6% for the He-r case and ~ 99% for the RG-r case), but only loses 18% of the mass. If the 
angular momentum, J = aM^^E^ui^, has a constant a after the impact of SN explosion, then a can 
be calculated using the angular momentum from the initial conditions (a = Jo/(M*i?^a;*) = 0.251). 
However, the post-impact radius of the companion is not straightforward to define, because the 
remnant companion star is not in hydrostatic equilibrium and the location of the photosphere is 
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unclear. We can estimate the equilibrium radius of the remnant star in the case MS-r using the 
virial theorem, -Rremnant ~ GM'^/Ei^t = 1-6 x 10^^ cm= 2ARq. Therefore, the post-impact rota- 
tional velocity can be calculated from the post-impact R' , M' and J' to obtain v'^ ~ 37 km s~^, 
which is only 23% of the initial rotational velocity (^^,0 = -R*w* = 164.4 km s~^). However, the 
virial theorem and constant a assume hydrostatic equilibrium, and Tycho's SN is only 439 yr old. 
The real R' should be larger, causing a lower rotational velocity for the post-impact companion 
star. This result provides additional support for the interpretation of Tycho G as the candidate for 
the progenitor star in Tycho's SN with the SDS. 

In contrast, the luminosity and surface temperature of the He star in the He- WD scenario are 
about 1.3 X 10^ Lq and 4.4 x 10^ K, which are too high to be consistent with the candidate progenitor 
for the Tycho's SN. This would also be the case for the RG-WD scenario since the luminosity of the 
remnant is expected to be comparable to its RG progenitor, which must be sufficiently luminous. 
This is required in the SDS model in order that the mass transfer rate prior to the SN explosion 
be greater than r^ 10^"^ MQyr~^ such that the CO WD accumulate sufficient matter to reach the 
Chandrasekhar limit. 



5. Conclusions 

We have investigated the impact of SN la ejecta on companion stars in the single-degenerate 
scenario via three-dimensional hydrodynamical simulations. We studied possible binary companion 
models, including a MS star, a RG and a He star, and considered the effects of asymmetry intro- 
duced by orbital motion and s pin. A detailed setup of the sub-grid SN la explosion using the W7 



model in lNomoto et al.l (| 19841 ) also is described. It is found that the orbital motion and spin lead 
to ~ 16% more unbound mass in the MS star companion channel but do not significantly affect 
the kick velocity. Furthermore, the orbital motion and spin play an important role in determining 
the morphology of the SNR. A power-law relation between the unbound mass and initial binary 
separation is found for all companion channels and is consistent with previous studies. Similarly, 
the kick velocity can be fitted by a power-law for the MS and He star binary companions. For the 
RG companion, we can only report a 40 km s~^ kick as an upper limit due to numerical uncertainty. 
By using the technique of passive particles, we find that the unbound mass is dominated by ablation 
instead of stripping. This result is in conffict with previous understanding, in which ablation was 
ignored in previous analytical studies. In addition, a hole shadowed by the ejecta is found to break 
the symmetry of the SNR. The amount of nickel contamination of the companion star is found to be 
~ W~^Mq for the MS star companion, ~ 10~^Mq for the RG companion, and ~ 10~^Mq for the 
He star companion. The corresponding nickel/iron to hydrogen plus helium abundance ratio may 
be useful for identifying the progenitor candidate in SN la remnants in future observations. We 
also find that the post-impact companion star loses about half of its initial angular momentum for 
the MS- WD scenario with the rotational velocity decreasing to 23% of its initial rotational velocity, 
providing added support for the SDS model for the Tycho SN. 
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The simulations presented here were carried out using the NSF Teragrid's Ranger system at 
the Texas Advanced Computing Center under ahocation TG-AST040034N. FLASH was developed 
largely by the DOE-supported ASC/Alliances Center for Astrophysical Thermonuclear Flashes at 
the University of Chicago. This work was partially supported by NSF AST-0703950 to Northwestern 
University and by the Computational Science and Engineering (CSE) fellowship at the University 
of Illinois at Urbana-Champ aign. Analysis an d visualization of simulation data were completed 



using the analysis toolkit yt (JTurk et al. 
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